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model 

R. StubnerE I. V. Tokatly, and O. Pankratov 

Lehrstuhl fiir Theoretische Festkorperphysik, Universitdt Erlangen-Numberg, Staudtstr. 7, 91058 Erlangen, Germany 

(Dated: 15 June 2004) 

We investigate the description of excitonic effects within time-dependent density-functional the- 
ory (TDDFT). The exchange-correlation kernel /xc introduced in TDDFT allows a clear separation 
of quasiparticle and excitonic effects. Using a diagrammatic representation for /xc, we express its 
excitonic part f^^ in terms of the effective vertex function A. The latter fulfills an integral equation 
which thereby establishes the exact correspondence between TDDFT and the standard many-body 
approach based on Bethe-Salpeter equation (BSE). The diagrammatic structure of the kernel in the 
equation for A suggests the possibility of strong cancellation effects. Should the cancellation take 
place, already the first-order approximation to fy^^ is sufficient. A potential advantage of TDDFT 
over the many-body BSE method is thus dependent on the efficiency of the above-quoted cancel- 
lation. We explicitly verify this for an analytically solvable two-dimensional two-band model. The 
calculations confirm that the low-order J^" perfectly describes the bound exciton as well as the 
excitonic effects in the continuous spectrum in a wide range of the electron-hole coupling strength. 

PACS numbers: 71.10.-w, 71.15.Q1, 71.35.-y, 78.20.Bh 
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I. INTRODUCTION 

Calculation of electronic excitation spectra remains 
one of the central problems of the quantum theory of 
solids. Of special interest are two-particle electron-hole 
excitations which determine the material's optical prop- 
erties. In semiconductors and insulators the electronic 
screening is suppressed by the energy gap and the in- 
teraction of the excited quasiparticles may substantially 
modify the excitation spectrum. The excitonic effects 
stemming from this interaction comprise the formation of 
bound electron-hole states as well as the alteration of the 
absorption in the continuum spectrum above the band 
edge. The latter is commonly referred to as unbound 
exciton effects or Sommerfeld absorption enhancement. 

In many-body perturbation theory two-particle excita- 
tions are characterized by the two-particle Green function 
which satisfies the Bethe-Salpeter equation (BSE)Ji Al- 
ready in 1980 Hanke and Sham^, using an approximate 
tight-binding representation, showed that the BSE cor- 
rectly describes the strong excitonic features above the 
optical absorption edge in Si. In the current "state of the 
art" procedure (see Ref. l3 for a recent review) the cal- 
culation of excitonic effects involves three steps. First, a 
density-functional theory (DFT) calculation in the local 
density approximation (LDA) is performed. On the sec- 
ond stage, the LDA Kohn-Sham (KS) energies and wave- 
functions are used as a starting point for the GW calcu- 
lation of the quasiparticle spectrum. Finally, the BSE 
is solved numerically, using the GW eigenvalues and the 
LDA wavefunctions as input characteristics of the nonin- 
teracting quasiparticles. The outlined procedure leads to 
highly accurate results, as has been shown for a number 
of relatively simple systems, mostly bulk semiconductors 
(see Ref. 3 and references therein). However, this method 
is extremely laborious, and for more complex systems the 
calculations become prohibitively expensive. 



A promising alternative, which is being intensively 
developed over recent years, relies on time-dependent 
density- functional theory (TDDFT) i4 This theory al- 
lows to calculate (formally exactly) the linear density- 
density response function and thereby the excitation 
energies,?: Since in the framework of DFT the exchange- 
correlation (xc) effects are lumped in a local xc poten- 
tial fxc, the TDDFT equation for the response function 
contains the variational derivative of v^c with respect to 
density /xc(r, i;r', t') = (5wxc(r, t)/(^n(r', t'). This xc ker- 
nel /xc is the central unknown quantity of TDDFT in 
the linear response regime. In their pioneering work 
Zangwill and Soven^ calculated the photo-absorption in 
rare gases in a self-consistent field manner. They used 
what later became known as adiabatic local density 
approximation (ALDA), simply substituting the time- 
dependent density in the LDA xc potential v^^^^ = 
v^^[n{r.t)). The resulting xc kernel is local in space 
and time f^^^^ = 6{v ~ v')5{t - t')dv^^^ / dn{r) . ALDA 
has been successfully applied to various finite systems 
like atoms or molecules i^^^iSiS Typically in these systems 
already the RPA response function calculated with KS 
eigenvalues and eigenfunctions gives good results. The 
correction due to f-^^^^ is quite small, which signifies 
that Hartree effects dominate in the response function. 
Unfortunately, f^tt^^ remains insignificant also in ex- 
tended systems like semiconductors or insulators, where 
KS-RPA gives a very poor description of the absorption 
spectraii2*ii Thus whereas a correct accounting for xc ef- 
fects becomes crucial in extended systems, the ALDA 
kernel /^(t^^ f^-i^^ to provide even a reasonable starting 
approximation. 

In the late nineties it has been realized that ALDA 
cannot serve as the basis approximation for the dynamic 
xc response of an inhomogeneous electron gas, because of 
the intrinsically nonlocal nature of fxc^ For extended 
systems with an energy gap this remains valid even in 



the static caseJ^^ii 

The importance of the nonlocahty of /xc was high- 
Ughted by the work of Reining and coworkersi^ti^ who 
were able to describe the contributions of unbound ex- 
citons in several diamond or zinc-blende type semicon- 
ductors with a static xc kernel proportional to l/|r — r'|. 
Other examples for this are the exact exchange kerne]ii£*i£ 
and the results of de Boeji et al<^, where in the context 
of time-dependent current-density functional theorjiSSiSi 
the nonlocal effects were crucial for accurately describing 
the effects of unbound excitonsiS^ 

It can be easily understood that a nonlocal /xc is cru- 
cial for describing excitonic effects. Within TDDFT the 
proper polarization operator is defined via the RPA-like 
equation 

X(w) = Xs{^) + Xs(w) • Uc{^) ■ xi^^), 

where in a crystalline solid x, the Kohn-Sham response 
function xs, and /xc are matrices in reciprocal space. The 
matrix structure of x is responsible for local-filed effects. 
However, these are relatively small in typical semicon- 
ductors and can be neglected for a qualitative analysis. 
Keeping only diagonal elements with zero reciprocal lat- 
tice vectors, we can easily solve for x obtaining 



x(^,q) = 



Xs(w,q) 



1 — lim 



l-/xc(tJ,q)xs('^,q)' 
The macroscopic dielectric function sm is given by 

euiuj) = 1 - lim Vc(q)x('^,q) 

Vc(q)xs(^,q) 

9~o l-/xc(w,q)xs(t^,q)' 

with the Coulomb interaction Vc(q) = Aiie^/q^. An 
additional excitonic peak in £m(w) appears when the 
denominator vanishes. However, it is well known that 
Xs is proportional to q^ in the limit g ^ for sys- 
tems with an energy gapiSSiSi Hence /xc must behave 
as l/q^ in this limit to counterbalance xs- Otherwise the 
xc kernel would have no effect on £m(w) at all. For the 
static long-ranged xc kernel of Reining et al. we have 
/xc((^, q) — Aire^P/q^ with some constant /3. The macro- 
scopic dielectric function thus reads 



Em — 1 — 



47re^as('-^) 
1 — 47re^/3as('^) ' 



where asC'^) is the macroscopic polarizability of the 
Kohn-Sham system as((^) = limg_,o Xs('^5q)/'Z^- For a 
typical as{uj) close to the band edge this formula sug- 
gests the existence of only one excitonic peak. However, 
one expects several peaks from unbound excitons above 
the band gap and bound excitons within the gap. Phe- 
nomenologically one could overcome this problem by in- 
troducing a frequency-dependent (3. One though would 
need to introduce very rapid oscillations in the region of 
the Rydberg series of bound excitonic states . 



Probably the most promising path in the quest for a 
good approximation to /xc is a direct comparison of the 
TDDFT formalism with the BSEJ^i2SiSL21 Simply com- 
paring the calculated spectra, it was found that it is often 
sufficient to use an approximation to /xc which is of the 
first-order in the screened particle-hole interaction. Al- 
though these results are very encouraging, it is unclear 
why this approximation is so efficient and what its range 
of validity is. 

In this paper we derive a diagrammatic expression 
which exactly relates the excitonic part of /xc to the BSE. 
We start with splitting /xc into two parts separately ac- 
counting for quasiparticle and excitonic effects. We then 
apply the diagrammatic rules we previously derived^^ to 
these two parts of /xc. This leads us to an expression for 
the excitonic part of /xc in terms of the three-point func- 
tion A. The latter satisfies an integral equation similar 
to the BSE which establishes the exact correspondence 
between TDDFT and common many-body theory. The 
main advantage of this approach is that the possibility 
of cancellation effects, which have been conjectured in 
Ref. 1271 is directly seen in the kernel of the equation for 
A. 

In order to investigate the properties of our integral 
equation and the applicability of low-order approxima- 
tions we study a model two-band system. In this model 
both the BSE and TDDFT equation can be solved an- 
alytically, which offers an ideal test bed for approxima- 
tions to the exact /xc. We find that indeed there are 
strong cancellation effects in the integral equation for A 
in the energy range close to the band gap. For this reason 
both the shallow excitons and the unbound excitonic ef- 
fects are well described with a first-order approximation 
to the excitonic part of /xc. 

The paper is organized as follows. In Sec. E] we inves- 
tigate diagrammatic properties of the excitonic part of 
/xc and derive an exact correspondence between TDDFT 
and BSE. In Sec. lIIII we introduce the model system that 
is used in further calculations. Sections IIVI and [V] are 
devoted to analytic calculations of excitonic effects with 
a general short-ranged and the Coulomb interaction re- 
spectively. Finally, we present our conclusions in Sec. lVIl 



II. DIAGRAMMATIC MEANING OF /xc 

The Bethe-Salpeter equation in the particle-hole chan- 
nel is commonly formulated as an integral equation for 
the particle-hole propagator or the scattering matrix T 
in the ladder approximation. This integral equation is 
equivalent to a summation of all ladder diagrams. The 
diagrammatic representation of the T-matrix formulation 
is depicted in Fig. Q^a) , where the full lines are the quasi- 
particle Green functions and the dashed lines represent 
the screened interaction. In the context of a connec- 
tion to TDDFT we are not interested in the four-point 
particle-hole propagator but rather in the response func- 
tions. Therefore we consider a modified BSE where two 
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FIG. 1: Diagrammatic representation of the Bethe-Salpeter 
equation 



(b) /E»= (^^ 

FIG. 2: (a) The n-th order of f^^ and (b)the integral equation 
for the three-point function A which is an essential part of f^^ 



of the external lines of the scattering matrix have been 
contracted to form the three-point function T. The dia- 
grammatic representation of the integral equation for T 
and its relation to the proper polarization operator x are 
depicted in Fig. ^b). For our purposes this equation is 
equivalent to the BSE. In the following we refer to the 
equation in Fig. ^h) as the Bethe-Salpeter equation. 

In terms of TDDFT the proper polarization operator 
X is given as 

X(w) = Xs{^) + Xs(w) • /xc(t^) • xi^^), (1) 

where xs represents the density-density response func- 
tion of the noninteracting KS particles, i.e., a bare loop 
of two KS Green functions. Equation ^ looks like 
the RPA equation, although it relates xs with the full 
X which includes all self-energy and ladder diagrams. 
We can therefore interpret /xc as an effective interaction 
which describes self-energy and ladder diagrams in the 
annihilation channel. Thereby /xc contains both quasi- 
particle and excitonic effects. As in this paper we are 
only interested in the excitonic effects, it is tempting to 
separate these two contributions to /xc, as suggested in 
previous works i-'^^i^^ This separation is indeed possible 
without approximations, because we have 



/x 



Xs' 



r' 



Xs'-Xqp+Xqp 



x-y 



(2) 
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where xqp is the density-density response function for 
the noninteracting quasiparticles. By definition, /^^ and 
f^^ are the kernels of the following RPA- type equations: 



XQp(w) = xs(w) + Xs(w) • f$^{uj) ■ XQp(w) (3a) 

X(w) = XQP('^) + XQP('^) • /x''c"('^) • X{^)- (3b) 

The newly introduced quantities f^^ and f^^ describe 
quasiparticle and excitonic effects respectively. This can 
be visualized by applying the diagrammatic rules for /xc 
as derived in Ref. |2^ The structure of the diagrammatic 
representation of /^^ and f^^ is similar to the one for 
/xc, except that for /^^ one has to use KS Green func- 
tions and should only account for all possible self-energy 
insertions in every order of the perturbation theory. This 



clearly describes the quasiparticle corrections. For f^^ 
one should use the quasiparticle Green functions with 
all possible particle-hole interactions. Obviously, this re- 
flects excitonic contributions. It is important to note, 
that the general properties of the perturbative expansion 
of /xc obtained in Ref . |2^ remain valid separately for f^^ 
and f^^. In particular, /^^ remains finite at KS excita- 
tion energies in every order of the perturbation theory. 
The same holds for f^^ at excitation energies of the non- 
interacting quasiparticle system. Note that our f^^ is 
the same as the /J^^ of Refs. 29 andM 

Let us briefly outline the diagrammatic rules for the ex- 
citonic part of /xc . According to Ref. |23 with the above- 
mentioned modifications we have to draw all loops with 
n particle-hole interactions to construct the n-th order 
correction. These diagrams serve as parent graphs for 
the construction of the n-th order f^^. To comply with 
the BSE we must use the ladder approximation here as 
well]^ Therefore, only one diagram with n interactions 
is left. To the two ends of the diagram we have to at- 
tach wiggled lines representing Xqp- Next, we work out 
all possibilities to separate this parent graph into two by 
cutting two fermionic lines. Then we join the external 
fermionic lines of these parts, connect them by a wig- 
gled line and change the sign of the resulting diagram. 
Obviously, the only way to separate the parent graph is 
to cut between adjacent interaction lines. The cutting 
does not change the ladder structure of the diagrams 
as seen in Fig. EI a). The summation of all ladder dia- 
grams can be cast in an integral equation as displayed 
in Fig. [2Ib). If one inserts f^^ obtained by solving the 
equation of Fig. EJb) into Eq. H3b|l and calculates the 
response function Xj the result will be the same as the 
X obtained from the BSE. In this sense the equation of 
Fig. |2Ib) gives the exact "translation" of the BSE into 
the TDDFT language. 

At this point we would like to highlight the differ- 
ences between the integral equation of Fig. I3b) and the 
equation for the excitonic part of /xc derived in Ref. |23 
(Eq. (4) and (5) therein). The iterative equation of 
Marini et al. is for the two-point xc kernel and is based 
on (finite order approximations to) the xc part of the 
response function and is logically analogous to our dia- 
grammatic expansion of /xc obtained in Ref. ^. How- 
ever, the equation of Fig. |2Ib) is an integral equation for 



(a) ^ -^^J>.v^v,x^ 

(b) T ^T -l5>^^^^ 



FIG. 3: Diagrammatic representation of the two possible in- 
terpretations of tiie difference between BSE and our equation 
for A. 



a three-point function analogous to the BSE^^ and can 
be solved instead of the BSE to obtain the same results. 

Unfortunately, the exact calculation of f^^ from the 
equation of Fig. |2Ib) is at least as difficult as obtaining 
an exact solution of the BSE. However, one can hope that 
the two-point kernel f^^ is more suitable for approxima- 
tions. An indication in this direction can be seen directly 
from the diagrammatic equation of Fig.[2Ib). Comparing 
this equation to Fig. ^h) we see that it can be obtained 
from the BSE by substitution of the particle-hole propa- 
gator as shown in Fig.PJa). This is the same replacement 
that was used in Ref. |28| to prove the cancellation of di- 
vergencies in /xc at KS excitation energies. Similarly, it 
facilitates the cancellation of divergencies in f^^ at QP 
excitation energies. 

Alternatively we can interpret the difference between 
the BSE and our equation for A as a modification of the 
interaction process in the second diagram of the r. h. s. of 
the equation in Fig.^b). This change of the interaction 
is displayed in Fig.|2Ib). As f^^ is responsible for trans- 
lating the ladder diagrams into the annihilation channel, 
it is nontrivial only if the ladder channel and the anni- 
hilation channel are distinguishable. This is easily seen 
from the replacement of Fig. Olb). Indeed the ladder 
and the annihilation channel coincide for nonrclativistic 
systems with a static point interaction. In this case the 
quasiparticle propagators in Fig.|2Ib) form a polarization 
loop which cancels the wiggled line. As a result the two 
diagrams of Fig. Ofb) cancel exactly. In the equation for 
A in Fig. ^h) this means that the last two terms cancel 
and A reduces to the first term. The excitonic part of 
the xc kernel then reduces to the interaction itself. One 
can therefore expect that in systems with a short-ranged 
and almost static effective interaction, the two terms of 
Fig-E^b) will cancel to a large extent, and a low-order ap- 
proximation to f^^ will be sufficient. Conditions like this 
can, e.g., be found in simple metals. On the contrary, in 
semiconductors screening is less efficient and the effective 
interaction is long-ranged. Further research is needed to 
verify to what extent the cancellation is efficient. 

Cancellation effects as we are expecting them from our 
diagrammatic equation have been seen in R ef. [ 271 the 
success of the lowest-order /xc found in Ref. |23 implies 
that de facto the cancellation can be efficient in materials 
with a band gap as well. 



III. MODEL SYSTEM 

In this section we consider a model system which re- 
veals both bound and unbound excitonic effects and 
where both the BSE and the equation for A can be solved 
analytically. With the exact f^^ at hand we can verify 
under what circumstances the first-order approximation 
to f^^ may be sufficient. The approximate kernel must 
describe bound as well as unbound excitons. A simple 
system with a bound exciton is given by the two-band 
Dirac model with a static density-density interaction. 
Moreover, we consider the two dimensional case in or- 
der to avoid technical difficulties with diverging integrals. 
The model Hamiltonian is given by 

+ - I d?r I d^r'n{Y)V{r - r')ri(r'), (4) 

where fi(r) — ■0'l'(r)V'(r) is the density operator, -(/'(r) is 
a two-component field operator, and V[r — r') describes 
the interaction between the particles. The Hamiltonian 
for the noninteracting particles reads 



A fc_ 
H = kxUx + kyUy + Actz = ( r a 



(5) 



± iky, the 



where (Jx,y,z are the Pauli matrices, k± = 
2D momentum operator k — (k^, ky), and the band gap 
equals to 2A. The energy dispersion for the noninter- 
acting particles has two branches which we label c and v 
for the (unoccupied) conduction band and the (occupied) 
valence band: 



£;c/v(k) = ±Ek = ±VA2 + fc2. 



The eigenvectors of H are 

*c. = ( .T 

with 



l-v 



Uk 



(6) 



(7) 




Note that this model can be understood as a two- 
component ( "relativistic" ) system. The interactions in 
the ladder channel and in the annihilation channel in such 
a system are always distinct regardless whether the inter- 
action is long- or short-ranged. Therefore in this model 
f^^ is nontrivial even in the case of a point interaction. 

We are going to solve this model in the ladder approx- 
imation, i.e., ignoring self-energy terms and higher-order 
corrections to the irreducible scattering matrix. Com- 
paring to the BSE this implies that H refers to the inde- 
pendent quasiparticles and V is the screened interaction 
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FIG. 5: The diagrammatic representation of the "bare" vertex 
and the "charges" of Eqs. 11211 and 11311 . 



FIG. 4: The real and the imaginary part of the exact xqp- 



between them. Therefore the one-particle Green func- 
tions of H are the quasiparticle Green functions: 



C?Qp(t^,k) 



*ck*L , *vk^lk 



u) — Ek + iS uj + Ek — iS 



(9) 



Note that ^E'vk and \I>ck are two-component vectors and 
Gqp is a 2 X 2 matrix. To lowest order in the wave vector 
q the quasiparticle response function xqpCi^iQ) is given 
by 

XQp('^,q) "" ^* / ^ / ■7^:^^'^C!Qp{e + u),k + q)GQp{e, 
0)2 -t- 1 , 1+OJ 1 



In- 



-5 - (10) 
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where uj = a'/(2A). The real and the imaginary part 
of this function are displayed in Fig. 01 A nonvanish- 
ing imaginary part occurs only at frequencies above the 
quasiparticle gap a) > 1, when the argument of the loga- 
rithm becomes negative. Note that one sees here explic- 
itly the g^-dependence of the response functions men- 
tioned in the introduction. 

In principle, we are now in the position to solve the 
BSE in this model analytically. However, to proceed fur- 
ther wc need to introduce some technical issues. Due to 
the matrix structure of Gqp one has to compute traces 
when calculating xqp. Similarly, the three-point func- 
tions r and A are 2x2 matrices and traces have to be cal- 
culated over internal indices in the diagrams of Figs.^a) 
and I3b). Since all these matrices have nonvanishing 
off-diagonal elements, evaluation of these traces becomes 
quite tedious. Therefore it is convenient to choose the 
eigenstates of H as the basis. The Green function Gqp 
becomes then a diagonal matrix with the elements 



c-c or v-v two-particle propagators vanish due to the in- 
tegration over frequency since both "upper" and "lower" 
Green function have their pole in the same half of the 
complex plane. This means that these diagrams do not 
contribute to the polarization. The same holds for the 
three-point functions F and A. Instead of one equation 
for the 2x2 matrix F as depicted in Fig. ^b) we obtain 
two coupled equations for the scalar functions Fcv and 
Fvc- For Fcv the upper line is a conduction-band state 
and the lower line a valence-band state, whereas for Fvc 
it is vice versa. 

The transfer to the diagonal representation of Gqp 
is equivalent to a corresponding transformation of the 
l.^ field operators in Eq. Q) . This transformation generates 
the "bare" vertex, describing interaction with an external 
field in the polarization diagrams, as well as the "inter- 
action vertices". In the original representation |(SJ all 
these vertices are simply unit matrices. They become, 
however, nontrivial matrices in the diagonal representa- 
tion. In fact, only three matrix elements of these vertices 
are essential. These are depicted in Fig. (Sj To the low- 
est order in the transfered wave vector q the bare vertex 
computes to 



Gc/v(w,k) = 



1 



Lu ^ Ek^iS 



(11) 



We thus can abandon the matrix notation altogether and 
use the two scalar Green functions Gc and Gv instead of 
the now diagonal 2x2 matrix Gqp. Now, every full line 
is either a conduction- or a valence-state propagator (|ll|l . 
This, of course, increases the number of diagrams we need 
to draw, because we have to consider all possible com- 
binations of conduction- and valence-band states. How- 
ever, in all diagrams the "upper" and "lower" Green func- 
tions that constitute a particle-hole propagator must al- 
ways be of a different type (c or v). The diagrams with 



7(q,k) 



*ck+q*vk 



1 

2Ek 



ulq- 



(«'¥)^ 



9+ 



(12) 
Note that 7(q, k) is linear in q because of the opposite 
parity of c- and v-eigenfunctions at k = 0. Thus the 
linear g-dependence stems from the off-diagonal momen- 
tum operator in the Hamiltonian (jSJ . This is actually the 
cause of the q^-dependence of the response function we 
refered to above in conjuction with Eq. H10|) and in the 
introduction. The "charges", i.e., the vertices associated 
with the interaction, are given by 

fc_ k' 

gi (k, k') = *];k*ck' = UkUk' + —VkVk' -^ (13a) 

k k' 

52 (k, k') = *ck*vk' = -^VkUk' - UkVk'-rr- (13b) 

All other possible combinations of valence- and 
conduction-band states differ from Eqs E| and El only 
by sign changes or complex conjugation. 

From now on we use a diagram technique with two 
different types of full lines representing conduction- and 
valence-band states. With these lines we associate the 
scalar Green functions of Eq. (|ll|l . Vertices are associ- 
ated with the scalar functions of Eqs. H12|) and \i'A\ . 

There is an alternative interpretation of this basis 
transformation. Gonsider one of the traces that has to 



be calculated for xqp 



tr(vI/ek+qVl/Tj^^^ 



)(*vk*:k), 



(14) 



where parenthesis show the grouping of the matrix mul- 
tiplication. One thus has to calculate the outer products 
of two vectors, multiply the resulting matrices and in 
the end take the trace. However, this grouping can be 
changed as follows 



tr*ck+q(*L+q*vk)* 



t 

vk 



*ck+c*vk)(*lk*ck+q) 



(15) 

Now one computes inner products of two vectors and 
multiplies the resulting scalar functions. Taking the 
trace is accounted for automatically in Eq. (|15|l . The 
change of the diagram technique outlined above is in ef- 
fect a way to incorporate this regrouping into the formal- 
ism. 

As noted above, working in the basis of the conduction- 
and valence-band states we have to split the BSE of 
Fig. nib) into two scalar equations for rcv(w,q, k) and 
rvc('-^, q, k). Note that while Fcv and Fvc depend on 
two momenta, they depend only on one frequency for 
a frequency-independent interaction. These two three- 
point functions are in fact not independent but related 
by the replacement q — > — q, w — > —lu and complex con- 
jugation. We can therefore derive one equation for Fcv: 



Fcv(c^,q,k)=Fi(w,q,k) (16a) 

.gi(k,k')Fcv(c^,q,k').9i*(k',k) 



j:y-:- 



2Ek' ~ uj 



^ ^^ g2(k, k')F,%(-c., -q, k')g2(k', k) 






M.,: 



Fi(L^,q,k)= V^k.k- 



2Ek' + uj 
.gi(k,k')7(q,k')gi*(k',k) 



k' 



/^ 



^^k,l 



2Ek' - uj 

g2(k,k07*(-q,kOg2(k^k) 
2Ek' + i^ 



(16b) 



Note that we omit the q-dependence in the "charges" gi 
and §2 as well as in the energy denominators, as we are 
only interested in the lowest-order expansion in q, which 
stems from the dipole matrix elements in the external 
vertices. Here 14;. k' are the matrix elements of the in- 
teraction between the particles. Now we investigate this 
equation for two different types of interaction, a general 
short-range interaction and the long-ranged Coulomb in- 
teraction. 



IV. A SHORT-RANGE INTERACTION 



The solution of the BSE 



scale shorter than A"-'^. The final results are expressed 
in terms of the physical (renormalized) scattering length, 
which includes the high energy contribution to the inte- 
grals in Eq. (|16|) . Having in mind this renormalization 
we can formally use a momentum independent bare in- 
teraction Vk.k' —V in Eq. H16|l . Note that this does not 
make f^^ trivially equal to the interaction itself as dis- 
cussed in the last paragraph of section^ because for our 
two-band model the annihilation channel and the ladder 
channel remain different even for a contact interaction. 
From Eq. p6|l we see, that the k-dependence of Fcv is 
given by the "charges" from Eq. (|13|l and has the same 
general form as for the bare vertex (|12() . The same is 
true for the q-dependence. Therefore we can write the 
following ansatz for Fcv^ 

Fev(c^,q,k) = ulq^T[^Xu) + {vk^)\+Ti^{uj). (17) 

Inserting this ansatz into Eq. 116|) we obtain two cou- 
pled equations: 



rff(")=V'i:^ 



2Ek V (2Sfc - uj) {2Ek + uj) 



k 



ul 



2Ek - UJ 



2Ek + ijJ 



r^t^. 



w 



r(^;*(--) 



(18a) 



^'S^-)-yY.^ 



2Ek V (2Sfc + w) (2Sfe - uj) 



yj: 



2Ek -UJ 



riiH^) 






(18b) 



Fcv (w) = Fcv (— w) and the equation to solve reduces 



From these equations it immediately follows that 
to 

1 f < 



v'tK^)^yY. 



2Ek V (2£^/c - uj) {2Ek + uj) 



yj: 



2Ek - UJ 2Ek + UJ 



ritH^ 



■.Vji{ij) + VKo{ij)Tl:l\uj). 



(19) 



At this point it is convenient to perform the above- 
mentioned renormalization of the interaction by splitting 
Kq into a low- and a high-energy part, the latter being 
equal to ^j. ^-g-. The high energy part logarithmically 
diverges at large k. This divergence can be removed by 
the standard renormalization of the interaction 



In this section we solve Eq. (I16|) for a short-ranged in- 
teraction, i.e., an interaction with a characteristic length 



V = 



V 



l-VE 



1 

k 4Ek 



47r 



(20) 





V C 

FIG. 6: The diagrammatic expression for IIxc 



Calculating the 2D integrals in Eq. H22() and 71 (w) in 
Eq. l(T^ we obtain 



47r 4cj 1 — uj 47r 



where we introduce the dimensionless scattering length 

a. With this renormalized interaction, Fct (w) fulfills the 
following equation 



r['^{Lu) = F7i(c^) + VKo{Lo)rg\Lo) (21) 



with 



i^oM-i^oH 



2^ Ar„_ Z^ ■ 



(2A + w) 



('221 



Equation H21|) has the obvious solution 



.is)., _ VliH 



nv{^) = 



l-VKo{uj) 



and 



71 (,w) 



167r V 2^2 ^ 1 _ jI, 


1 


16. .i(^^") 2)- 





(25) 



Inserting these integrals into the solution (|21() we can 
compute the response function's xc part IIxc — X ~ Xqp 
and its first-order approximation IIxc ■ In diagrammatic 
(23) form IIxc is displayed in Fig. El The results are 



and 



nxc('i,q) = - 



ni'c^(^,q) 



167rAi 



(FiCo) 



\f 



167rAw2 I i_aF(w) ^ l-aF*{~Cj) 



{{F{u)-lY + {F*{^u)-lY) 



(26) 



(27) 



r 



where the function F{ui) is defined in Eq. H24|l . 

The equation for A given in Fig. |2tb) can be solved in 
a similar fashion, which allows then to calculate the exci- 
tonic part of the exact xc kernel. However, it is easier to 
obtain f^^ directly from Eq. |(2J) and the exact response 
fmiction x = Xqp + IIxc with IIxc from Eq. (|26|) : 

Ex.- . XQp(^.q)nxc(^,q)XQM^,q) 

i + XQp(^>q)nxc(t^,q) 

From here the first-order approximation to f^^ immedi- 
ately follows: 

/x^,^^'^(^, q) = XQp('i, q)nW(i, q)xQ^(ti, q). (29) 

/•Ex(l) 
/xc 

imate solution for the response function's xc part 

niV(c:',q) 



Inserting this f^^ in Eq. (|3b|) we arrive at an approx 
on fc 

n^c"('^,q) = 



t(i)^,^, 



(30) 



i-n^i^((i',q)xQp((:',q) 

Note that although f^^ is based on IIxc , this formula 
does not coincide with Ilic . In accordance with Eq. (|3b|) 



it accounts for an infinite series of diagrams instead. This 
way, the excitonic pole in IIxc, Eq- l|26|l . which has been 
lost in nlc\ Eq. |2Z|), reappears in Eq. (jSOJ. 



B. Results 

Having calculated the exact and the approximate ex- 
pressions for the excitonic part of the xc kernel and for 
the response function we are now in the position to com- 
pare these results. Let us start with the excitonic peak 
in the absorption spectrum. In both the exact and the 
approximate response functions the excitonic peak origi- 
nates from the divergence of the xc part. The exact Hxc 
in Eq. (|26|l has a pole at 



1 - aF{Cj) = 



(31) 



for positive a and < d) < 1. When the dimensionless 
exciton binding energy e — \ — lj \s small, Eq. I|31|) has 
an approximate solution 



2 exp ( — 
a 



(32) 
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FIG. 7: The real and the imaginary parts of the exact x (full) 
and the x based on f^^ (dashes) for the scattering length 
a — 0.25. The real and the imaginary part of xqp are shown 
by the dotted line. 
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FIG. 8: The real and the imaginary parts of the exact x (full) 
and the x based on f^^ (dashes) for the scattering length 
a = 0.6 



The approximate response function of Eq. (|30|l which is 
based on f^^ has a pole 

l-XQp('^,q)n«(ci,q) = 0. 

Similarly to the solution (|32|l of the exact equation (|31|l , 
for a small binding energy e this equation can be approx- 
imately solved by 



e = 2 exp — 
\ a 



1 1 - a - Vl - 2a - a^ 



2a 



(34) 



which differs from the exact result H32|) . However, e.g. 
at a = 0.2 this is only about 14 % larger than the exact 
solution 132() . which gives £ « 0.013. For a smaller scat- 
tering length and therefore a smaller binding energy the 
agreement between the exact result and the one based on 
f^^ is even better. 

In Fig. 13 the real and the imaginary part of the exact 

X and the x based on f^^ are shown for the scatter- 
ing length a = 0.25. For comparison, xqp is also dis- 
played. One clearly observes a very good agreement be- 
tween the exact and the approximate response functions. 
Both main features, the enhancement of the imaginary 
part (the Sommerfeld factor) and the excitonic peak are 
correctly reproduced. The latter is indicated by the ar- 
rows on the plots of the imaginary part. The apparent 
difference between the exact and the approximated re- 
sponse function is the exact position of the excitonic peak 
as described above. The oscillator strengths should also 
be different, which is not, however, reflected in the figure. 
From the above discussion one could conjecture that 
the approximation for the response function based on 



FIG. 9: The real and the imaginary part of the exact f^^/a 
for a — 0.05 (dots), a = 0.1 (dashes), a — 0.25 (long dashes), 

and a — 0.6 (dot-dashes) compared with f^^ /a (full). 



f^^ is always sufficient. This is, of course, not true, 
as we can see from Fig. |S1 where the approximate and 
the exact response function are compared for the larger 
value of the scattering length a = 0.6. With this large 
interaction^ the position of the excitonic peak in the 
approximate response function is clearly wrong. In addi- 
tion, the imaginary part has the wrong magnitude above 
the band gap. 

To uncover the background of the good agreement be- 
tween the exact and the approximate response functions 
for a "weak" interaction, let us compare f^^ with the 
exact f^^ for different interaction strengths. We note 



first that f^^ is proportional to the scattering length 

a. Therefore, more adequate is to compare f^^ /a with 

(1) 



f33) f^cl^i 3'S done in Fig. |H| It is clearly visible that f^^ 

is a very good approximation to the exact f^^ in a fre- 
quency range close to the band gap. Actually, for uj = 1 
we obtain from our exact expressions 



f^f\^ = l)--'-^a 



and 



/xc ('i - 1) 



167rA 



l-a2/4' 



(35) 



(36) 



Thus even for a very strong interaction of a = 0.6 



where /, 



Ex(l) 



leads to the rather poor response function 



Ex(l) 



at o) = 1. In 



of Fig. |S1 there is only a 10 % error in /, 

the static case the errors are larger, as we get for a) = 0: 

(37) 



.Ex(i).~ n^ 167rA9 
/xc (w = 0) = ^-a 



q 



and 



/xc ('^ = 0) 



167rA9 



1 



(38) 



Here a 10% error is already reached for a — 0.1, as the 
ratio of f^^ and f^^ is of order a, whereas it is of order 
a^ &t Cj — 1. 

It is also interesting to look at the difference 5f^^ be- 



tween the exact f^^ and /, 



pEx(l) 



strengths. Since the leading order term in 5f^^ is of 



for different interaction 

Ex 



^Sf^^la- 
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the order a^, we normalize these differences by a^, when 
plotting them in Fig. ^| The general behavior of the 
different curves is quite similar, which indicates that it is 
mostly the second-order approximation f^^ which con- 
tributes to Sf^^. One may be surprised that the curves 
for stronger interaction are closer to for w < 1, which 

might indicate that f^^ is a better approximation for 
stronger interactions. However, this is not the case. This 
only tells us that f^^ is positive for uj < 1. 



C. The validity of the first-order approximation 



From the previous section we can conclude that f^^ 
is a good approximation to f^^ close to the band gap 
(D = 1 practically for any interaction strength. Even for a 
strong interaction, where the response function stemming 
from f^^ is quite wrong, the first-order kernel f^^ is 
still very good close to di — 1. This is the reason for the 

success of f^^ in describing the bound excitonic states 
as shown in the previous section. Since around uj = 1 

the first-order kernel f^^ is a good approximation, if 
the (exact) bound excitonic state lies within this region. 



f^^ will describe it correctly. If, however, the bind- 
ing energy is outside this region, f^^ will fail. Note 
that as seen from Fig. this region gets smaller as the 
interaction increases and at the same time the binding 
energy of the exciton increases. Hence the error in the 
exciton binding energy increases with the increase of the 
interaction strength. 

Can we understand the good agreement between f^^ 
and f^^ for d) « 1 in terms of the integral equation for A? 
For this we explicitly calculate the diagrams involved in 
the replacement shown in Fig.|3Jb). Working with sepa- 
rate Green functions for the valence- and the conduction- 
band states we have four possible combinations of the ex- 
ternal lines and can therefore split this replacement into 
four parts according to these combinations. At energies 
close to the band gap the two-particle propagator with a 
conduction-band state in the upper line should dominate. 
Hence, the most important part of this replacement is the 
diagram with a conduction-band state in the upper line 
on both sides of the graphs. This diagram together with 
its translation into quantities introduced in the previous 
section is displayed in Fig. ^2 Note that we use the bare 
interaction V here, as the interaction renormalization in 
these diagrams means simply a replacement V ^ V. 



ck-|-q ck'-l-q ck-|-q 

-r 
I _ 

I 

V k V k' V k 






c k'-l-q 
vk' 



= VgfiK k') - ri(<.., q, k)xQU'^, q)7*(q: k') 

FIG. 11: Exphcit calculation of one replacement. Here to and 
q are the transfered energy and momentum. 

Inserting gi, Fi, xqPj ^^^id 7 in the expression in Fig. llll 
and taking the limit a) ^ 1 we obtain: 



^^^(k^k') _ri(cD = i,q,k)xQM'^ - i,q)7*(q,k') = 

A \ /A:_ \2/A:^ 



/ 9 9 / A \ /A:_ \2/^+ \^ fc- k' „/fc' \2/9-\2 A \ 



(39) 



r 



Not surprisingly the whole expression (|39|l is proportional 
to the interaction V (or V after renormalization) . In the 
proportionality factor in round brackets all summands 
contain (1 — A/Ek') (or powers thereof). For the first 
summand this is directly visible, in the others it is "hid- 
den" in Vk' ■ Let us look at this factor more closely. When 
the diagram of Fig. ^] is inserted into the equation for A 
of Fig. Efb) the integration over k' has to be performed. 
The main contribution to this integral comes from small 
momenta, due to the small energy denominators in the 



particle-hole propagator. However, for these small mo- 
menta the expression (I39f) is small as A w Ek> ■ 

We can therefore conclude that close to the band gap 
the kernel in the equation for A is indeed small for those 
states which mainly contribute to the integral. Thus we 
explicitly observe the cancellation we qualitatively dis- 
cussed in sectional This cancellation explains the excel- 
lent agreement between f^^ and f^^ for a) ss 1. 

A complimentary interpretation can be obtained from 
looking at the diagrammatic expansion of the response 
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(a) <Q) (b) (T><T><y> 

FIG. 12: Third-order term in the expansion of (a) the exact 
X and (b) the x based on f^^ in terms of the interaction. 



cc cv vc vv cc 



vv vc cv cc vv 
FIG. 13: Approximation for IIxc in the limit of small i. 



functions. In Fig. 1121 the third-order term in the expan- 
sion of the exact x and the x based on f^^ is displayed. 
The difference between the two expressions is similar to 
the replacements discussed in section ^ This remains 
true in all orders of the perturbation theory. From the 
previous calculation of the effect of this replacement fol- 
lows, that the exact % and the x based on f^^ are 
almost identical in the vicinity of the band gap, indepen- 
dently of the position of the excitonic peak. 



in which we are interested in here. The Is wavefunction 
needed in this approximation is given by 



V'o(?-) 



exp{-y/Aeor} 



and the exciton binding energy is 



(44) 



(45) 



Introducing dimensionless variables e and £o — £o/(2A), 
we arrive at 



V. THE COULOMB INTERACTION 

In this section, we consider the 2D Dirac model of 
Eq. (0J) with the Coulomb interaction between the parti- 
cles, i.e.. 



V{v- 



r — r' 



(40) 



where e is the particle's charge. Note that Eq. 140|) is 
the 3D interaction although our model system is 2D. We 
will not attempt to solve exactly the BSE with this in- 
teraction, but rather focus on the properties of shallow 
excitons. The weak binding limit e = 1 — a) <C 1 is, in 
fact, a "nonrelativistic" limit where the BSE reduces to 
the two-particle Schrodinger equation.'^S- The quasiparti- 
cle energy eigenvalues are approximately 



Ek-A + 



and for the eigenvectors holds 



2A' 



Uk 



1 and Vk ~ 0. 



(41) 



(42) 



Solving the BSE thus reduces to solving the positronium 
problem in 2D for particles with a mass A, i.e., for a 
reduced mass A/2. The response function can be written 
in the spectral representation as^i 



\^n{r =_0)|2 

k - 

oj-{2A~ eo) ' 



X(u;,q) = |7(q,0)p^ — . 

^^^ UJ - (2A - En) 



l7(q,o)| 



(43) 



where the 'i/'n are the eigenf unctions with eigenvalues 
En (positive for bound states) of the above-mentioned 
Schrodinger equation. The approximation in Eq. (|43f) is 
valid close to the excitonic peak of the Is ground state. 



x(e, q) 



£o 



27rA eo - e 



(46) 



for e close to the excitonic peak of the Is ground state. 
We now want to compare the position and the oscilla- 
tor strength of the excitonic peak of the Is ground state 



with the results obtained from f^^ . 
expand xqp of Eq. ifTUI) for small e 

,. . ci" (. (2 

XQp(e,q) 



For this we first 



167rAV^"U'"' 



(47) 



In order to obtain /, 

t(i) 



Ex(l) 



we also need an expression 
for IIxc'' in the limit of small e. In diagrammatic form 
IIxc is given by a single graph as displayed in Fig. ^1 
There are four possibilities to distribute the conduction- 
and valence-band states in the two particle-hole propaga- 
tors. The particle-hole propagators with a conduction- 
band state in the upper line and a valence-band state 
in the lower line diverge at ci) = 1. The particle-hole 
propagators with a valence-band state in the upper line 
and a conduction-band state in the lower line diverge 
at cj ~ —1. Therefore, in the limit where the energy 
goes to the band gap (or equivalently e is small), the 

main contribution to IIxc comes from the diagram with 
both particle-hole propagators of the upper line belong- 
ing to the conduction-band state and the lower line to 
the valence-band state. Neglecting the other contribu- 
tions is our first approximation in calculating IIxc . The 
other approximation is using Eqs. \M\ and H42() . With 
these approximations we obtain 



n«(£,q) 



l7(q,0)P 
(2A)^ 



■E 

k.k' 



F(k-k') 



fe\^ 



Hi) 



> A 



^ +£ 



(48) 
where F(k— k') = 27re^/|k— k'| is the Fourier representa- 
tion of V(v — r') from Eq. H40I) . The easiest way to solve 
this double integral is to look at it in the real space, in- 
stead of the Fourier space, where it becomes a single inte- 
gral. The 2D Fourier transformation of the particle-hole 
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propagator (i(fc/A)^ + e) gives the modified Bessel 
function of the second kind Kq , so that we can write 



niV(e,q) 



87rA\/2l 



dpK^p). 



(49) 



The integral over Kq computes to (7r/2)^, and we arrive 



at 



n«(£~,q) 



q^ TT^ 
IGttA 2 



(50) 



fEx(l) 



Now we have all ingredients to build /^ 

XQpHxc Xqp ^^"i insert it in Eq. Ij3bp . This gives the 
following approximation for the response function based 



onfl^^'^: 



x^ (£, q) 



XQp(£,q) 



i-niV(e,q)xQp(£~,q) 



(in(i: 



167rA 



(ln( 



(51) 



Note that this equation is valid only for small e, since we 
derived xqp and IIxc only for small values of e. 

From Eq. (|51|) we clearly see that X'^ contains an 
additional pole where the denominator vanishes. This is 
the excitonic peak in this approximation and its energy 
Eq is defined by 



(52) 



To evaluate the approximation used let us reformulate 
Eq. ^ as 





(53) 



We see that the approximate binding energy s'q and the 
exact binding energy Eq are identical for 



£0 = 2 exp 



1 



1 
189' 



(54) 



In other words, for the interaction strength which corre- 
sponds to this quite realistic exciton binding energy, the 
approximation using f^^ gives the exact binding en- 
ergy. For stronger and weaker interaction strength there 
is some error. To be more precise, the error is below 10 % 
for £o between 1/165 and 1/222. It is below 20% for Eq 
between 1/148 and 1/271. Note that in the range of en- 
ergies where the approximation gives good results, e is 
indeed small, so that our approximation is consistent. It 
is interesting to note that for the Coulomb interaction — 
in contrast to a short-range interaction — one does not 

get the correct exciton binding energy from f^^ in the 



limit of the interaction strength going to zero. Note also 
that while the exact solution of the 2D hydrogen prob- 
lem gives rise to an infinite series of excitonic peaks in 
the exact response function, we obtain only one excitonic 
peak in the approximate approach. This problem was al- 
ready touched in the introduction. Although f^^ is 
frequency dependent, it does not contain the rapid oscil- 
lations needed to describe the whole series of excitonic 
states. 

To compare the approximate response function of 
Eq. (|51|l to the exact one from Eq. (|46|l and determine 
the residual it is best to expand x around Eq. This can 
be done by calculating the first-order Taylor expansion 
of the denominator in Eq. (|51|l : 



,f(i) 

x^ (e, q) 



,2 -^^ ir 



IGttA 1 

2 



(55) 



In the case where the exciton energy is exactly repro- 
duced, i.e. for Eg = £o, the second fraction in Eq. H55|) 
is about 16.6. The oscillator strength is thus too large 
by a factor of two inspite of the fact that the binding en- 
ergy is exact. The value of this second fraction is almost 
constant for a large range of ratios £o/£o, so that mainly 
the Eq in the numerator of the third fraction in Eq. I|55|) 
produces additional errors in the oscillator strength. The 
larger oscillator strength is not too surprising, as the ex- 
citonic oscillator strength is distributed over more peaks 
in the exact response function. 

Similar to the results for the short-range interaction 
of the previous section, we observe that the cancellation 
effects are also effective for the long-ranged Coulomb in- 
teraction. However, the general structure of this cancella- 
tion is different and cannot be universally characterized. 

Hence the accuracy of f^^ has to be checked for every 
particular system. The integral equation for A together 
with the replacement procedure of Fig. |2| provides an ap- 
propriate tool for this task. 



VI. CONCLUSIONS 

In conclusion, we have investigated the excitonic effects 
on the response function within TDDFT. We have split 
the xc kernel into a quasiparticle f^ and an excitonic 
part f^^. Using a diagrammatic expansion we derived 
the integral equation for the three-point function A, in 
terms of which f^^ can be exactly expressed. As this 
integral equation is similar in structure to the BSE, it 
establishes the connection between the common many- 
body theory and TDDFT. The kernel of the equation for 
A shows the possibility of cancellation effects. If these 
cancellations were complete, f^^ would be equal to its 



first-order approximation /, 

fEx(l) 



Ex(l) 



This suggests that in 



some situations f^^ provides a good substitute to /, 



Ex 
xc ■ 
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We have presented explicit calculations for a model 
two-band semiconductor with short-ranged interaction 
which can be solved analytically. Comparing the ex- 
act response function with the response function derived 

from f^^ we confirm a very good agreement for a weak 
interaction, where both the position of the shallow exci- 
ton and the Sommerfeld factor are correctly described. 
We were able to trace this to the strong cancellation oc- 
curing for energies close to the band gap (almost) inde- 
pendently of the interaction strength. Calculations with 



the Coulomb interaction give similar results, though the 
deficiencies of f^^ are somewhat worse in this case. 

These calculations represent an example that the inte- 
gral equation for A can serve as a tool for evaluation of 
the validity of low-order approximations to f^^. Testing 
the predictive power of this approach for other systems 
will be the subject of further work. 

This work was supported by the Deutsche Forschungs- 
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